#include <vector>	//vector
#include <limits>	//min, max
#include <algorithm>	//sort
#include <sstream>	//stringstream
#include <ostream>	//ostream
#include <istream>	//istream
#include <string>	//getline
#include <map>		//map
#include <iostream>	//cin
#include <utility>	//pair

#include <cmath>	//abs, sqrt
#include <cstdlib>	//rand, srand
#include <ctime>	//time

#include "RPTree.h"
#include "RPLeaf.h"

using namespace std;

// Vide
RPTree::RPTree()
: v(0), usedRule(RPTree::none), threshold(numeric_limits<double>::signaling_NaN()), mean(0.0), depth(0), left(NULL), right(NULL)
{}

// Depuis une sauvegarde
RPTree* RPTree::load( istream &in )
{
	RPTree* root = NULL;
	string buffer;
	istringstream iss;

	map<void*, RPTree*> id;
	RPTree* curentRPTree;
	void* curentId;
	double curentData;
	int curentRule;

	while( in.good() )
	{
		getline(in, buffer);
		iss.clear();
		iss.str(buffer);

		iss>>curentId;
		iss.ignore(); // virgule

		getline(iss, buffer, ',');
		if( buffer.compare("NODE") == 0 )
		{
			id[curentId] = curentRPTree = new RPTree();
			if( ! root )
			{
				root = curentRPTree;
			}

			iss.ignore(); // [
			iss.ignore(); // virgule

			while( iss.peek() != ']' )
			{
				iss>>curentData;
				curentRPTree->v.push_back(curentData);
				iss.ignore(); // virgule
			}

			iss.ignore(); // ]
			iss.ignore(); // virgule
			iss>>curentRule;
			curentRPTree->usedRule = (RPTree::rule) curentRule;

			iss.ignore(); // virgule
			iss>>curentRPTree->threshold;

			iss.ignore(); // virgule
			iss>>curentRPTree->mean;

			iss.ignore(); // virgule
			iss>>curentRPTree->depth;

			iss.ignore(); // virgule
			iss>>curentId;
			curentRPTree->left = (RPTree*) curentId;

			iss.ignore(); // virgule
			iss>>curentId;
			curentRPTree->right = (RPTree*) curentId;
		}
		else
		if( buffer.compare("LEAF") == 0 )
		{
			unsigned int idLeaf;
			unsigned int depth;
			vector<double> vData(0);

			iss>>idLeaf;

			iss.ignore(); // virgule
			iss>>depth;

			iss.ignore(); // virgule
			iss.ignore(); // [
			iss.ignore(); // virgule

			while( iss.peek() != ']' )
			{
				iss>>curentData;
				vData.push_back(curentData);
				iss.ignore(); // virgule
			}
			
			id[curentId] = new RPLeaf( depth, idLeaf, vData);
		}
	}

	for( map<void*, RPTree*>::iterator i = id.begin() ; i != id.end() ; i++ )
	{
		if( i->second->usedRule != RPTree::none )
		// Pas une feuille
		{
			i->second->left = id[i->second->left];
			i->second->right = id[i->second->right];
		}
	}
	return root;
}

// Racine
RPTree::RPTree( RPTreeParam &param, list< list<unsigned int> > &leavesContent )
: v(param.ncol()), usedRule(RPTree::none), threshold(numeric_limits<double>::signaling_NaN()), mean(0.0), depth(0), left(NULL), right(NULL)
{
	vector<unsigned int> dataIdx(param.nrow());
	for( unsigned int i = 0 ; i < param.nrow() ; i++ )
	{
		dataIdx[i] = i;
	}

	unsigned int label = 0;
	srand ( time(NULL) );
	init( param, dataIdx, leavesContent, label );
}

// Noeud
RPTree::RPTree(const RPTreeParam &param, vector<unsigned int> &dataIdx, list< list<unsigned int> > &leavesContent, unsigned int depth, unsigned int &id)
: v(param.ncol()), usedRule(RPTree::none), threshold(numeric_limits<double>::signaling_NaN()), mean(0.0), depth(depth), left(NULL), right(NULL)
{
	init( param, dataIdx, leavesContent, id );
}

// Initialiseur
void RPTree::init(const RPTreeParam &param, vector<unsigned int> &dataIdx, list< list<unsigned int> > &leavesContent, unsigned int &id )
{
	vector<unsigned int>::size_type data_num = dataIdx.size();	// Nombre de données
	unsigned int data_col = param.ncol();				// Nombre de colonnes

	{
		//on tire un vecteur aléatoire et on le normalise
		double norme = 0.0;
		for( vector<double>::size_type i = 0 ; i<data_col ; i++ )
		{
			v[i] = ( ( ( (double) rand() ) / ( (double) RAND_MAX ) ) * 2.0 ) - 1.0;
			norme += v[i]*v[i];
		}

		norme = sqrt(norme);

		for( vector<double>::size_type i = 0 ; i<data_col ; i++ )
		{
			v[i] /= norme;
		}
	}

	//projection
	double delta;
	vector< pair<double, unsigned int> > a(data_num);

	{
		double min = numeric_limits<double>::max();
		double max = numeric_limits<double>::min();
		for( vector<unsigned int>::size_type i = 0 ; i<data_num ; i++ )
		{
			a[i].second = dataIdx[i];

			for( vector<double>::size_type j = 0 ; j<data_col ; j++ )
			{
				a[i].first += v[j] * param.get( dataIdx[i] ,j );
			}

			min = min > a[i].first ? a[i].first : min ;
			max = max < a[i].first ? a[i].first : max ;
			mean += a[i].first;
		}

		delta = max - min;
		mean /= data_num;
	}

	double deltaA= 0.0;
	for( vector<unsigned int>::size_type i = 0 ; i<data_num ; i++ )
	{
		deltaA += ( a[i].first - mean ) * (a[i].first - mean );
	}

	deltaA *= 2.0 / data_num;

	//Rule
	vector<bool> rule(data_num, false);
	unsigned int accTrue = 0;

	if( delta * delta <= param.c * deltaA )
	{
		vector<double> aAcc(data_num);
		vector<double> aAcc2(data_num);

		usedRule = RPTree::hyperplane;
		sort(a.begin(), a.end() );

		for( vector<unsigned int>::size_type i = 0 ; i<data_num ; i++ )
		{
			aAcc[i] = a[i].first + ( i ? aAcc[i - 1] : 0 );
			aAcc2[i] = a[i].first * a[i].first + ( i ? aAcc2[i - 1] : 0 );
		}

		vector<unsigned int>::size_type imin = 0;
		double min = numeric_limits<double>::max();
		for(  vector<unsigned int>::size_type i = 0 ; i<data_num-1 ; i++ )
		{
			double mu1 = 1.0 / (i+1) * aAcc[i];
			double mu2 = 1.0 / ( data_num - (i+1) ) * ( aAcc[data_num - 1] - aAcc[i] );
			double ci  = ( aAcc2[i] - (i+1) * mu1 * mu1 ) + ( (aAcc2[data_num - 1] - aAcc2[i]) - ( data_num - (i+1) ) * mu2 * mu2 );

			if( ci < min )
			{
				min = ci;
				imin = i;
			}
		}

		threshold = ( a[imin].first + a[imin+1].first ) / 2.0;

		for( vector<unsigned int>::size_type i = 0 ; i<data_num ; i++ )
		{
			accTrue += ( rule[i] = a[i].first <= threshold );
		}
	}
	else
	{
		usedRule = RPTree::distance_from_the_mean;
		vector<double> dist_from_mean(data_num);

		for( vector<unsigned int>::size_type i = 0 ; i<data_num ; i++ )
		{
			dist_from_mean[i] = abs( a[i].first - mean );
		}

		sort( dist_from_mean.begin(), dist_from_mean.end() );

		threshold = dist_from_mean.size()%2 ?	  dist_from_mean[ dist_from_mean.size()/2 ] :
							( dist_from_mean[ dist_from_mean.size()/2 - 1 ] + dist_from_mean[ dist_from_mean.size()/2 ] ) / 2;

		for( vector<unsigned int>::size_type i = 0 ; i<data_num ; i++ )
		{
			accTrue += ( rule[i] = abs( a[i].first - mean ) <= threshold );
		}
	}

	// Fils		
	vector<unsigned int> leftDataIdx;
	leftDataIdx.reserve(accTrue);
	vector<unsigned int> rightDataIdx;
	rightDataIdx.reserve(data_num - accTrue);

	for( vector<unsigned int>::size_type i = 0 ; i<data_num ; i++ )
	{
		
		if( rule[i] == true )
		{
			leftDataIdx.push_back(a[i].second);
		}
		else
		{
			rightDataIdx.push_back(a[i].second);
		}
	}

	if( leftDataIdx.size() != data_num and rightDataIdx.size() != data_num )
	// Les branches sont équilibrées
	{
		if( leftDataIdx.size() <= param.min_size or
		    depth+1 == param.maxdepth )
		{
			left = new RPLeaf( param, leftDataIdx, leavesContent, depth+1, id );
		}
		else
		{
			left = new RPTree( param, leftDataIdx, leavesContent, this->depth+1, id );
		}

		if( rightDataIdx.size() <= param.min_size or
		    depth+1 == param.maxdepth )
		{
			right = new RPLeaf( param, rightDataIdx, leavesContent, depth+1, id );
		}
		else
		{
			right = new RPTree( param, rightDataIdx, leavesContent, depth+1, id );
		}
	}
}

RPTree::~RPTree()
{
	if( left ) { delete left; }
	if( right) { delete right; }
}

RPTree::operator string() const
{
	ostringstream out;
	out<<this<<",NODE,[,";

	for( vector<double>::const_iterator i = v.begin() ; i != v.end() ; i++ )
	{
		out<<*i<<',';
	}

	out<<"],"<<usedRule<<','<<threshold<<','<<mean<<','<<depth<<','<<left<<','<<right;
	return out.str();
}

void RPTree::save( ostream &out ) const
{
	out<<(string) *this<<endl;
	left->save(out);
	right->save(out);
}

unsigned int RPTree::use( double *data ) const
{
	double scalProd = 0.0;

	for( unsigned int i = 0 ; i < v.size() ; i++ )
	{
		scalProd += v[i] * data[i] ;
	}

	switch( usedRule )
	{
		case RPTree::hyperplane:
			if( scalProd <= threshold )
			{
				return left->use(data); 
			}
			else
			{
				return right->use(data); 
			}

			break;

		case RPTree::distance_from_the_mean:
			if( abs( scalProd - mean ) <= threshold )
			{
				return left->use(data); 
			}
			else
			{
				return right->use(data); 
			}

			break;
		default:
			// impossible
			return (unsigned int) -1;
	}
}

void RPTree::use(const CVGreyc::NDArray::CArray<double>& ndArray, unsigned int *result) const
{
	for( unsigned int i = 0 ; i<ndArray.shape(0) ; i++ )
	{
		result[i] = use( ndArray.data() + ndArray.index()(i,0) ) ;
	}
}

const RPLeaf *RPTree::get_leaf_by_id(unsigned int id) const
{
	const RPLeaf *leafRef = left->get_leaf_by_id(id);
	return ( leafRef ? leafRef : right->get_leaf_by_id(id) );
}
